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Summary 

Consider the problem of estimating the median of N items to a precision e, i.e. the estimate (0. should be such that, 
with a large probability, the number of items with values smaller than \y is less than — ( 1 + |e| ) and those with values 



N 

greater than \i is also less than - (1 + |e|) . Any classical algorithm to do this will need at least Q. 
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Quantum mechanical systems can simultaneously carry out multiple computations due to their wave like properties. 
This paper gives an ^ j~j J ste P algorithm for the above problem. 

1 Introduction 

1.1 Statistical Problems Statistical problems involve processing large amounts of data. By their nature quantum 
mechanical systems can simultaneously carry out multiple computations and if adjusted properly, the overall compu- 
tation can be significantly faster. [Shor94] used the fact that the periodicity of a function can be obtained rapidly by a 
quantum mechanical algorithm to develop a logarithmic time integer factorization algorithm. [Gro96] presented an 

O ( JN) step quantum mechanical algorithm to examine an unsorted database containing N items, for a single item 
that satisfies a particular property. The next step is to explore what statistics can be obtained by quantum mechanical 
algorithms, faster than is possible with classical algorithms. 

Consider an experiment that can give two results with slightly different probabilities of Q + e J and 

i — e J respectively where the value of e is unknown. The experiment is to be repeated several times in order to esti- 
mate e. If the experiment is repeated v times, it can be easily shown that, with a large probability, the number of 
occurrences of the first event will lie in the range v (^ + e j± O (Vv) . Therefore in order for the uncertainty due to 
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Therefore in order to estimate the median with a precision of e (i.e. the number of items with values smaller than |i is 

' 1 ' 



O(Jv) to become less than v|e| , O (Jv) must be smaller than v|e| , equivalently v must be greater than Q. 



less than — (1 + |e|) and those with values greater than |i is also less than — (1 + |e|) ) will need Q. 



samples 



[Wilks43]. 

1.2 Quantum mechanical algorithms A good starting point to think of quantum mechanical algorithms is 
probabilistic algorithms [BV93] (e.g. simulated annealing). In these algorithms, instead of the system being in a spec- 
ified state, it is in a distribution over various states with a certain probability of being in each state. At each step, there 
is a certain probability of making a transition from one state to another. The evolution of the system is obtained by 
premultiplying this probability vector (that describes the distribution of probabilities over various states) by a state 
transition matrix. Knowing the initial distribution and the state transition matrix, it is possible in principle to calculate 
the distribution at any instant in time. 

Just like classical probabilistic algorithms, quantum mechanical algorithms work with a probability distribu- 
tion over various states. However, unlike classical systems, the probability vector does not completely describe the 
system. In order to completely describe the system we need the amplitude in each state which is a complex number. 



The evolution of the system is obtained by premultiplying this amplitude vector (that describes the distribution of 
amplitudes over various states) by a transition matrix, the entries of which are complex in general. The probabilities 
in any state are given by the square of the absolute values of the amplitude in that state. It can be shown that in order 
to conserve probabilities, the state transition matrix has to be unitary [BV93]. 

The machinery of quantum mechanical algorithms is illustrated by discussing the three operations that are 
needed in the algorithm of this paper. The first is the creation of a configuration in which the amplitude of the system 

being in any of the 2" basic states of the system is equal; the second is the Fourier transformation operation and the 
third the selective rotation of different states. 

A basic operation in quantum computing is that of a "fair coin flip" performed on a single bit whose states 



1 1 
1 -1 



. A bit in the state is 



are and 1 [Simon94]. This operation is represented by the following matrix: M = — 
transformed into a superposition in the two states: [ — — , ]. Similarly a bit in the state 1 is transformed into 

{ — \z\, i.e. the magnitude of the amplitude in each state is but the phase of the amplitude in the state 1 is 

inverted. The phase does not have an analog in classical probabilistic algorithms. It comes about in quantum mechan- 
ics since the amplitudes are in general complex. This results in interference of different possibilities as in wave 
mechanics and is what distinguishes quantum mechanical systems from classical systems. 

In a system in which the states are described by n bits (it has 2 n possible states) we can perform the transfor- 
mation M on each bit independently in sequence thus changing the state of the system. The state transition matrix rep- 
resenting this operation will be of dimension 2"x2". In case the initial configuration was the configuration with all 

n 

n bits in the first state, the resultant configuration will have an identical amplitude of 2 in each of the 2" states. This 

is a way of creating a distribution with the same amplitude in all 2 n states. 

Next consider the case when the starting state is another one of the 2" states, i.e. a state described by an n bit 
binary string with some Os and some Is. The result of performing the transformation M on each bit will be a superpo- 
sition of states described by all possible n bit binary strings with amplitude of each state having a magnitude equal to 



2 1 
2 and sign either + or -. To deduce the sign, observe that from the definition of the matrix M, i.e. M = — 

J2 



1 1 
1 -1 



the phase of the resulting configuration is changed when a bit that was previously a 1 remains a 1 after the transfor- 
mation is performed. Hence if x be the n-bit binary string describing the starting state and y the n-bit binary string 

describing the resulting string, the sign of the amplitude of y is determined by the parity of the bitwise dot product of 

x and y , i.e. (-1) . This transformation is referred to as the Fourier transformation [DJ92]. This operation is 
one of the things that makes quantum mechanical algorithms more powerful than classical algorithms and forms the 
basis for most significant quantum mechanical algorithms. 

The third transformation that we will need is the selective rotation of the phase of the amplitude in certain 



states. The transformation describing this for a 3 state system is of the form: 



e 
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e 











where i = J-i and 



e 

(j)p ^3 ^ arbitrary real numbers. Note that, unlike the Fourier transformation and other state transition matrices, 

the probability in each state stays the same since the square of the absolute value of the amplitude in each state stays 
the same. 



2(a) The Abstracted Problem Consider the problem described in the following figure: there are a number of 
states S 1 ,S2,-S N , each of which has a value V(S) associated with it. 

V(S) < H A V(S) > H 



N/2(l+e) T (i N/2(l-e) 

Figure Given a threshold u, estimate e, i.e. the difference in the number of states with values above & below [I. 

It should be clear that if we are able to solve the problem described in the above figure, then by repeating it with 
different u, a logarithmic number of times, we can obtain the median. Next consider the problem described in the fol- 
lowing paragraph in which a range for e is specified and it is required to estimate e in this range with a specified pre- 
cision 0. An algorithm for this problem will be presented in this paper. A little thought will show that the problem in 
the above figure, and hence the problem of finding the median, can be easily transformed to this problem (this is dis- 
cussed further in the appendix). 

Let a system contain N = 2" states which are labelled S 1 ,S 2 ,-S N . These 2 n states are represented as n bit strings. 
Let each state, S, have a value, V(S), associated with it. A quantity (X is given that partitions the items such that the 

number of states, S, with value, V(S), less than u, is equal to — (1 + e) , where |e| is known to be in the range 

O.IEq < |e| < £q for a specified (where Eq is assumed to be smaller than 0.1). The problem is to estimate e with a 

precision of (where is positive) i.e. if E { be the estimate of e , then with a probability approaching 1, the true e 

1 



lies in the range: e 1 (l-0)<e<e 1 (l + 0) . The algorithm of this paper solves the problem in O 



e o 2 



time 



steps. The problem of finding the median is easily transformed to this, as indicated in the appendix (incidentally, even 

i 

to solve this abstracted problem classically, will need a time Ql —— 

e 2 / 



for the same reason indicated in section 1.1). 



2(b) Some notation In addition to the symbols defined in section 2(a) above, the following notation will be used: 

(i) In all steps of the algorithm, the amplitudes of all states, S, with values V(S) < (J. are the same, this amplitude will 
be denoted by k; similarly in all the steps of the algorithm, the amplitudes of all states, S, with V(S) > (J. are the same, 
this will be denoted by /. The amplitudes in states, S, with V(S) < (J. will be referred to as the "k amplitudes"; the 
amplitudes in states, S, with V(S) > (X will be referred to as the "I amplitudes". 

(ii) The heart of the algorithm is the loop in step l(iv) of the algorithm in the next section. This increases the k ampli- 
tudes in each repetition of the loop. This is done in a way so that the increase is approximately proportional to e. After 
V repetitions of the loop, the value of the k amplitude is denoted by fe y (e) & the value of the / amplitude is denoted 

by / v (e) . 

(iii) As mentioned in section 1.2, the probabilities are given by the square of the absolute values of the amplitudes. 
Hence the sum of the squares of the absolute values of the amplitudes in all states must be unity. However, in this 
paper we take the sum of the squares of the absolute values in the states to be N, where N is the number of states - nor- 
malized amplitudes can be obtained by dividing the amplitudes given in this paper by . 

(iv) O(n) & £2(n) (this is standard notation) A function f(n) is said to be O(ri) if there exists a constant K such that 
for all n greater than some specified n Q , f(n) < Kn . Similarly a function f(n) is said to be (n) if there exists a 

constant K such that for all n greater than some specified n„ , f{n) > Kn . 



3 Algorithm This algorithm solves the problem described in 2(a): given (a, e Q , , where e ft , > , it finds an 



-()' 



estimate of e such that if e L be this estimate, then with a high probability, the true e lies in the range: 
EjO-0) <e<e 1 (l + 0) . 



(1) Repeat the steps, (i)....(iv), a = O 



(i) 



(ii) 



(iii) 



(iv) 







times. 



Initialize the system to the distribution 

( 1 , 1 , 1 .... 1 ) , i.e. there is the same 
amplitude to be in each of the N states. This 
distribution can be obtained in O (log AO steps 
as discussed in section 1.2. 
In case the system is in a state with V(S) > |i, 

rotate the phase by ^ radians. 

Apply the shift transform S which is defined by 
the NxN matrix S as follows: 

S = TV + TV > 11 Q > S = TV - A — TV - 

pq N N PP N V N 

This shift transform, S, can be implemented as 
S = FRF , where F is the Fourier Transform 
Matrix as discussed in section 1.2 and R is: 

R = if p^q&R = 1 if p = 
pq pp ' 

& R = -i if p * . 
pp 1 

As discussed in section 1 .2, the Fourier 

transformation matrix F is defined as 



F = 2- n/2 

pq 



(_1) P ■ 9 ) where p is the 



binary representation of p , and p ■ q is the 

bitwise dot product of the n bit strings p & q . 
Repeat the following steps - (a), (b), (c) & (d), 

" - ii 

(a) 



'0 



times. 



In case the system is in a state, S, with 
V(S) < |i, rotate the phase by % radians, 
(b) Apply the diffusion transform D which is 
defined by the matrix D as follows: 



D 



pq 



- if p * q & D 
N PP 



■4 



This diffusion transform, D, can be 
implemented as D = FTF , where F is 
the Fourier Transform Matrix, as 
discussed in (iii), and Tis defined as: 

T =0 if p±q & T = 1 if p = 

pq pp 1 



& T 



PP 



-1 if/7^0. 



A A 



Step (l)(i): States with values above & below p., 
have the same amplitudes. 



A A 



>' > 

Step(l)(ii): States with values above (J, have the 

phases of their amplitudes rotated by 90 
The amplitude vector becomes ( 1, 1, i, i). 



A n 



Step(l)(iii): The system shifts to states with 
values above (0.. The amplitude 
vector becomes (0, 0, 1+i, 1+i). 

(Note that in case the number of states on the 

two sides is different, then there will be some 

amplitudes in states with V(S) < 



In each repetition of loop (l)(iv)(a)...(d) 
the k amplitude increases by £2(e). 

Therefore in P repetitions, the k amplitude 
becomes £2(Pe). 



(c) In case the system is in a state S, with 
V(S) > |A, rotate the phase by % radians. 

(d) Apply the diffusion transform D as 
described in (b) above. 

(v) Sample and record the resulting state. 

(2) The relative probability of the system being in a state 
If, ,2 

with V(S) < n is -(J fep (e) | ( 1 + e) J . In case the 

fraction,/ of the a samples (obtained in step l(v) 
above) be calculated, then with a confidence level 

greater than ( 1 - 2exp^-2K 2 j j , the quantity 
^ |fep (e) | 2 ( 1 + e) j will lie in the range (f- K0) 
to (f + K0) (theorem 8). Knowing these bounds for 



The k amplitude (which is £2(Pe)) is 
estimated by taking a samples (a is 
defined in step 1). 



the function jfcp (e) | ( 1 + e) ; e can be bounded 
f K0A 

within a range of -g- J with the same confidence 

level (theorem 9). As defined in step l(iv), 
1 



P = Qy-J ■ 11 follows that the maximum expected 
error in e will be O (K0e) . 



4. HOW does the algorithm WOrk? Step l(i) initializes the system so that there is an equal amplitude to 

n 

be present in all of the N = 2 states. In step l(ii), states with value greater than u, have the phases of their ampli- 

71 

tudes rotated by - radians. Next, in step l(iii), a shift transform is carried out by which the system shifts to the states 

with value greater than |i, i.e. the k amplitudes are approximately zero and the / amplitudes are of magnitude approx- 
imately J2 . In case there were an equal number of states with values above and below u, the amplitude in states with 
value smaller than \y will be zero. In general, if the difference in the number of states above and below \y is Ne, the 
amplitudes in states with values smaller than p. will be proportional to e. This is somewhat interesting since such a 
redistribution could not be carried out classically. However the probabilities are proportional to the squares of the 
amplitudes; therefore if we try to estimate e by sampling the distribution and counting the number of samples with 

f n 

values below u, it will take £1 — samples before we obtain even a single sample with value below u. As discussed 

in the introduction (section 1.1), this time is of the same order as would be required to estimate e classically. 

This paper presents a technique by which the magnitudes of the k amplitudes can be increased by an additive 

constant of £2(|e|) in each iteration sequence. Therefore by repeating this iteration sequence Cl( — ^] times, the 



magnitudes of the k amplitudes can be made to be Q — . By repeating this process Q. 



si „m 



times and observing 



the fraction of samples with values less than u, — can be estimated with a maximum error of O (0) and thus |e| can 

£ 



be estimated with a maximum error of (e Q 0) . The heart of the algorithm is this iteration sequence for increasing 

the magnitudes of the k amplitudes by Gl (|e|) and by a quantity proportional to |e| . This process is carried out in 
steps l(iv)(a)...(d) of the following algorithm, its logic is discussed in the rest of this section. 

One key step in the loop in steps l(iv)(a)...(d) is the diffusion transform D. By means of theorems 4(a) & 4(b) 
and by using superposition, it is possible to study the effect of the diffusion transform on arbitrary values of k & I. 
Theorem 4 (a) - In case the initial amplitudes be k = 1,1 = 0, then after applying the diffusion transform D, the 
amplitudes become k = £ , I = \ + £ . 

Theorem 4 (b) - In case the initial amplitudes be k = , / = 1 , then after applying the diffusion transform D, the 
amplitudes become k = 1 - e , I = -£ . 

Using theorems 4(a) & (b), it is possible to study the evolution of arbitrary initial amplitudes by using superposi- 
tion. Let the initial distribution before steps l(iv)(a)...(d) be given by k = K & I = L . 

(a) After step l(iv)(a), the distribution becomes k = -K & I = L. 

(b) The distribution after step l(iv)(b) is obtained by theorems 4(a) & (b) and by superposition. Applying this 
gives: k = - Ke + L (1 - e) & / = - (1 + e) K- eL . 

(c) The distribution after step l(iv)(c) is obtained by changing the signs of the / amplitudes. This gives: 
k = -K£ + L(l-£) &/= (1 +e)K+eL. 

(d) The distribution after step l(iv)(d) is obtained by theorems 4(a) & (b) and by superposition. Applying this 

gives: k = A"f 1 - 2e 2 J + Z^2e - 1£ J & / = - K^2e + 2e 2 ] + A 1 - 2e 2 J . This is theorem 7. 

Theorem 7 - In case the initial amplitudes before step l(iv)(a) of the algorithm be k & / , then after steps 



l(iv)(a)...(d), the new amplitudes, k & I are k = K\ 1 - 2e )+L\2e-2e J and I = Ky- 2e - 2e )+L\l-2e 
respectively. 

As in [BBHT96], the evolution of the amplitudes in successive repetitions of the loop can be followed by the follow- 
ing corollary: 

Corollary 7.2 - For the algorithm of section 2, the k amplitudes after r repetitions of the loop in step l(iv) are given 

2 2e-2e 2 

by k r = 7(1 + e + /) sinr<|) + ecosr(|) , where 7 = and (]) is defined by the following equation: 

2e + 2e 

2 

cos* = 1 - 2e . 



By definition, cos(|) = 1 - 2e ; therefore for small e , (]) = sin(f> = 2e . For small values of r(|) , sinr(|) = r(|) = 2r£ & 
cosr(|)= 1 . 7, as defined in cor. 7.2 above, is approximately equal to 1 for small 8. Therefore by cor. 7.2, 

\k\ = 2 V2r|e| . Therefore in r = Cl( - J repetitions of the loop of l(iv), \k\ = £1 ( 1) and it is also proportional to 



|e| . Thus by repeating the experiment described in step 1 of the algorithm Q. 
with a maximum relative error of O ( 0) . 



( 1 1 



times, it is possible to estimate |e| 



5 Proofs The following section proves that the system discussed in section 3 is indeed a valid quantum mechani- 
cal system and it gives the desired result with a probability O(l). The proofs are divided into three parts. First, it is 
proved that the system is a valid quantum mechanical system; second, that in each repetition of the loop of step l(iv), 
the magnitudes of the k amplitudes increases by £2 (|e|) ; finally it is proved that by calculating the fraction of the 
J 1 1 



a = O 



samples with values less than |i , the parameter |e| can be estimated with a precision of . 



5(a) Preliminary results In order to prove that the algorithm of section 3 is a valid quantum mechanical system 
we need to prove that D & S are unitary and can be implemented as a product of local transition matrices. The diffu- 
sion transform D & the shift transform S, are defined by the NxN matrices D & S, as follows: 

(5.0) D = J, if p*q&D = -1 + 2 

PI N PP N 

(5-1) S = i + l ifp^q & S = I-j[^=i 

PI N N PP N V N 

Since D & S can be implemented as products of unitary matrices (as required in steps l(iv)(b) and l(iii) respectively, 

and as proved in theorems 2(a) and 2(b)), both D & S can be seen to be unitary. However, since a direct proof is short 

and simple, a few lines will be devoted to this. The proof for unitarity of the matrix D is the same as in [Gro96]. 

Theorem 1(a) - D is unitary. 

Proof - For a matrix to be unitary it is necessary and sufficient to prove that the columns are orthonormal. The magni- 

n 

tude of each column vector, q, is ^ , q = 1,2. ..n. Substituting from (5.0), it follows that this is 

p=l 

n 

2^2 4 

1 + (N - 1) — which is equal to 1. The dot product of two column vectors is > D. D. . This is 

NJ N 2 4 F L-i jp jq 

7=1 

4 4 

N 2 N 

Theorem 1(b) - S is unitary. 

n 

Proof -The magnitude of each column vector, q, is ^ ^ q > 1 = l,2...n. Substituting from (5.1), it follows that this 



p=l 



2 

is — ^ + — }^ + 2 — r~ — which is equal to 1. The dot product of two distinct column vectors is S S 



where S denotes the complex conjugate of S . Evaluating > S S gives: 
rq r J ° rq ° ' ■ rp rq ° 

r = 1 

1 .(JV-l)Yl i\ fl iYl . (AT- 1) ^ (AT- 2) .... 

i— — — — + - + -- + 1 — tt-'- + 2- 1 which is equal to zero. 

N N AN NJ In NAN n ' n 



The way D & S are presented above, they are not local transition matrices since there are transitions from each state to 
all N states. Using the Fourier transformation matrix as defined in section 3, they can be implemented as a product of 
three unitary transformations, each of which is a local transition matrix. 

Theorem 2(a) - D can be expressed as D = FTF , where F is the Fourier Transform Matrix as mentioned in section 

1.2, and/? is defined as follows: T = if p*q &T = 1 if p = & T = -1 if p±0. 

pq 1 "* pp 1 pp 1 

Proof - We evaluate FTF and show that it is equal to D. As discussed in section 3, F = 2~ n/2 (-1) P ' 1 , where p 
is the binary representation of p , and p ■ q denotes the bitwise dot product of the two n bit strings p and q . T can be 
written as T = + T^ where 7^ = -/ , / is the identity matrix and ^ 00 = ^ ' T 2 pq = if /? ^ 0, g ^ . It is 
easily proved (it follows since the matrix M, as defined in section 1.2, is its own inverse) that FF=I and hence 

= FT^F = -/.We next evaluate D^ = FT^F . By standard matrix multiplication a( j = ^.^ntJo bc^cd' 

be 

Using the definition of T 2 , it follows that D n , = 2F n F „ , = — (-1) ~ a ' + ' ~ d = - . Thus all elements of the 

& z ' 2, ad aO Oil 2» N 

2 

matrix D 2 equal — . Adding Dj and D 2 the result follows. 



Theorem 2(b) - S can be expressed as S = FRF , where F is the Fourier Transform Matrix as mentioned in section 
1.2, and R is defined as follows: R =0 if p^q&R = 1 if p = & R = -i if p * . 

pq pp pp ' 

Proof - We evaluate FRF and show that it is equal to S. As discussed in section 3, F = lr n ^ (-l)P ' a , where p 
is the binary representation of p , and p ■ q denotes the bitwise dot product of the two n bit strings p and q . R can be 
written as R = R^+ R^ where R^ = -il , / is the identity matrix (with ones on the diagonal) and ^ 2 00 = ^ + ' ' ' 

/? 2 = if p^0,q^0 .As mentioned in theorem 2(a), FF=I and hence Sj = FR^F = -il . We next evaluate 

S 2 = FR^F . By standard matrix multiplication ad = ^^ a ^2 bc^cd ' Using me definition of R 2 , it follows 

fee- 
that S, , = (1 +i)F „F„ , = (1 +;) (_i)5-0 + 0-a = (! + Xhus all elements of the matrix 5 2 equal 

2, ad v 7 flO Od 2 n N 

— ^-^ . Adding 5; and 5 2 the result follows. 



5 (b) Loop of Step l(iv) Theorems 5,6 & 7 prove that the magnitudes of the amplitudes in states with values 
smaller than (0. (i.e. the k amplitudes) increase by 0(|e|) as a result of each repetition of the loop in steps 

l(iv)(a)...(d). Hence in repetitions of the loop, the amplitudes of those states will attain a value 0(1). 

Theorem 3 - After step l(iii) of the algorithm, the state vector is as follows: k = e , / = ( 1 + e) + i . 

Proof - After step l(i), the amplitudes in all states are equal, i.e. k = 1,1 = 1 ; after step l(ii), the / amplitudes have 

their phase rotated by - radians and the amplitudes become k = 1 , / = i . The matrix S for the shift transform is an 

NxN square matrix defined as follows: S = - + -, if pl=q&S = - - i\ ^-— T - 
H pq N N PP N \ N 

The amplitudes k & I after step l(iii) thus become: 

MM¥)M?<^-0(^K<'-^H 
<-<M¥MS«-'>-0(^KM^ -»♦•>♦' 



Theorem 4 (a) - In case the initial amplitudes are k = 1,1 = 0, then after applying the diffusion transform D, the 
amplitudes become k = £ , / = 1 + E . 

Proof -As mentioned before (4.0), the diffusion transform D is defined by the JVxJV matrix D as follows: 
(5.0) D =-, if p±q &D =-! + -. 

pq N pp N 

Therefore the amplitudes after applying D become: 

*-('-M?»«>-0-« 

I = f'l-f ' (1 + e) = d + e). 

Theorem 4 (b) - In case the initial amplitudes are k = , I = 1 , then after applying the diffusion transform D, the 
amplitudes become k = 1 - e , / = -e . 

Proof -As in theorem 4(a) above, the amplitudes after D become: 
2 N 

k = W 1 T (1 " e) = (1_e) 

'-(■-!M?<'->-0 — 



As is well known, in a unitary transformation the total probability is conserved - this is proved for the particular case 
of the diffusion transformation. This enables us to follow the evolution of the algorithm while keeping track of only 
one of the two quantities k or I, not both. 

Theorem 5 - In case the initial amplitudes be k = K, I = L, then after applying the diffusion transform D, the new 

2 2 2 2 

amplitudes k = K^, I = satisfy ( 1 + e) + Lj ( 1 - e) = K ( 1 + e) + L ( 1 - e) . 

Proof -Using theorems 4(a) & 4(b), it follows by superposition that the new amplitudes after the diffusion transform, 

(2 2 \ 

D, are given by K [ = Ke + L ( 1 - e) & = K ( 1 + e) -Le . Calculating \K { ( 1 + e) + Lj ( 1 - e) I , the result 
follows. 



In order to follow the evolution of the k & I amplitudes due to successive repetitions of the loop in step l(iv) (which 
is done in theorem 7), we first prove theorem 6: 

Theorem 6 (a) - In case the initial amplitudes before step (a) of the algorithm be k = 1,1 = 0, then after steps (a), 

2 2 

(b), (c) & (d), the amplitudes are k = 1 - 2e , / = -2e - 2e . 

Proof - After step (a), the amplitudes become fe = -l,/ = 0.It follows by theorem 4(a) that after step (b) the 
amplitudes become fc = -e,/ = -(l+e); after step (c) the amplitudes become k = -E, I = ( 1 + e) . Using the- 
orem 4(a) & 4(b), it follows by superposition that after step (d), the amplitudes are 

k = -e- e+ (1 + e) • (1-e) = l-2e 2 & I = -e(l +e)-(l + e)e = -2e-2e 2 . 



Theorem 6 (b) - In case the initial amplitudes before step (a) of the algorithm be k = , / = 1 , then after steps (a), 

2 2 

(b), (c) & (d), the amplitudes are k = 2e - 2e , I = 1 - 2e . 

Proof - After step (a), the amplitudes stay unchanged as k = , / = 1 . It follows by theorem 4(b) that after step 
(b) the amplitudes become k = 1 - e , / = -e ; after step (c) the amplitudes become k = ( 1 - e) , / = E . Using 

2" 



theorems 4(a) & 4(b), it follows that after step (d), the amplitudes are k = E • ( 1 - e) + £ • ( 1 - e) = { 2e - 2e ] & 
I = ( 1 + E) ( 1 - E) - E • E = 1 - 2e 2 . 



Finally using theorem 6, we can follow the evolution of the k & I amplitudes due to successive repetitions of the loop 
in step l(iv). 

Theorem 7 - In case the initial amplitudes before step (a) of the algorithm be K & L , then after steps (a), (b), (c) & 

(d), the new amplitudes, k & I are k = Al- 2£ 2 J + l(^2e - 2£ 2 j and I = A- 2e - 2z J + zi 1 - 2£ 2 
Proof - Using theorem 6(a) & 6(b), it follows by superposition that the amplitudes after steps (a), (b), (c) & (d) are 

given by k = K 1 - 2£ 2 + L 2e - 2£ 2 & I = - K 2e + 2£ 2 + L 1 - 2£ 2 



As in [BBHT96], the evolution of the state vector is obtained in terms of the sin & cos functions: 



Corollary 7.1 - In case the initial amplitude be k = ysinx and / = cost where y = 



2 2e-2e 



; then after r repe- 



2e + 2e 



titions of the loop l(iv), k^ = ysin + x) , / = cos (r(|) + x) where cos(|) = 1 - 2e . 

Proof - We prove the result by induction. By the initial conditions, the result clearly holds for r = . Assume it to be 



true for a particular r . We show that it holds for (r + 1) . By theorem 7: 



V+ 1 



7+1 



1-2e 2e-2e 

2 2 
-2e-2e 1-2e 





k 




r 




I 




r 



By 



the induction hypothesis, k r = ysin (r<|> + x) , / = cos (r<\> + x) . Expressing the transformation matrix in terms of 
y and § from the definitions in the statement of the corollary, it follows that: 



r+ 1 



r+ 1 



cos(j) ysin(|) 
--sind) cosd) 

y 



ysin (r(|) + 1) 
cos (r(|) + x) 



. By trigonometric identities for the sum of the sine and cosine functions: 



cos(|) ysin(|) 

--sind) cosd) 
Y 



ysin (r<|) + x) 
cos (r(|) + x) 



ysin ((r+l)<|) + x) 
cos ( (r+ l)(|) + x) 



Hence 



>+ 1 



V+ 1 



ysin ((r+l)(()+x) 
cos ( (r+ l)(|) + x) 



Thus by 



induction, the corollary holds for all r. 



Corollary 7.2 - For the algorithm of section 2, the k amplitudes after r repetitions of the loop in step l(iv) are given 
by k r = y( 1 + e + i) sinr(|) + ecosr(|) . 

Proof - The initial amplitudes are k = 8 & / = 1 + E + ; . The result follows by cor. 7.1 and superposition. 

1 lei 
Corollary 7.3 - There exists a number P < j^- such that |fep| > 



'40e, 







20e n 



Proof - By definition (in cor. 7. 1), cos(|) = 1 - 2e ; therefore for small £,(]) = sin<f> = 2e . Also for small values of 
r(|) , sin r§ = r§ & cos r(|) = 1 ; assuming small e , it follows from the definition of y in cor. 7.1 that y = 1 . Therefore 

by cor. 7.2, \k \=2j2r\e\. 



5(c) Estimating 8 : In case the probability of an event is estimated by carrying out a identical experiments, 
it follows by the law of large numbers, that the maximum expected error in the estimate is o( — ] . By the argument 

yJaJ 

before cor. 7.3, the amplitude of the system being found in states, S, with V(S) < is Q. ((3e) , the probability is pro- 



portional to the square of the amplitude and is hence £il P e j = Q. 



,2 2 



( 2\ 
8 



since the experiment of section 2 is 



repeated Q. 



( 1 ) 

ve j 



times, the maximum expected error in [ — ] is (9(0) ; i.e. K-O(B) <( — ] <K+0(Q) 



where K = 0(1) (this is proved in theorem 8). Equivalently E (K { -O (0)) < |e| < £q (K^ + (9(0)) where 
= (9(1), hence e can be estimated with a precision of (9(0) as desired in section 2(a) (proved in theorem 9). 



In order to prove theorem 8, we first state a version of the well known law of large numbers from probability theory: 
Lemma 8 - In an experiment, let the probability of the event A be p, the probability of the complementary event A 



will be (1 —p) . In T| repetitions of the experiment, the fraction of occurrences of A lies between p-K^J- and 

P + kJ^ with a probability of at least 1 - 2exp^-2K 2 
Proof - Omitted. 



Theorem 8 - In case a = O 



' — 



samples be obtained (as described in step 1 of the algorithm) and a fraction /have 



a value less than u, then with a confidence level greater than 1 - 2exp^-2K 2 j j , the quantity ^ \k^ (e) | 2 ( 1 + e) 
lies in the range (f- K0) to (f+ K0) . 



1*0 (e) I 

Proof - By the definition of (e) , it follows that the probability of each state with V (S) < |i is H ^ . Since 

there are ^ ( 1 + e) such states, it follows that the combined probability of all these states is ^ jfcp (e) | 2 ( 1 + e) 

Using lemma 8, with the event A being the event corresponding to the occurrence of a sample with value smaller than 
H , the theorem follows. 



Next we show that the value of e can be derived with a maximum expected error of O (08) from the estimate of 
i ,2 

|fcp (e) | ( 1 + e) , which was obtained in step 2 of the algorithm from the fraction of samples with values less than 
(0. . For this we prove a preliminary lemma: 

Lemma 9 - Let a and b be two arbitrary real numbers such that a<b . Let Af = f(b)- f(a) and the derivative/'^ 

at every x, such that a < x < b , is greater than d (where d > ). Then b -a<^ . 

d 

Proof - By the remainder theorem, f(b) -f(a) = f (x) (b - a) , for some x such that a<x<b . Therefore 
f{b) -f{a) >d(b-a) . 



Theorem 9 - Let the function /(e) be ^ (e) | 2 ( 1 + e) J . Let \f(b) -f(a) | = 2k0 where P« < 0.1 , $b < 0.1 . 

Then \b-a\ < 10, OOOk0£ q . 

Proof - Follows by these steps: 

2 2 2 2 

(i) By cor. 7.2, |fep (e) | = |y ( 1 + e) sinfity + ecosP<|)| + y | sin P<f»| ,(|) = 2e, ysl-e. 

(ii) ^(|^ p (e)| 2 (l +e)) = A[2(sinf3^) 2 J = A[ 8 (pe) 2 J = 16p 2 e. 



(iii) By lemma 9, \b-a\< ^ (b) ~ f{a) ^ . Therefore \b - a\ < = = 

min\f(x)\ 8P e 4P 8 4P e Q 



f 2^ 
2 



e< ,(10) e. 



ve J 4 — 
20 



1 
(the last step in (iv) follows by using cor. 7.2 by which P > ^q^ - & by tne problem specification by which — < 10 ) 

Therefore in case the value of e be estimated from the fraction of samples obtained in step l(iv) with values less than 

2 1 



\i (as described in step 2); then with a confidence level greater than \ 1 - 2exp^ -2k J J , e will lie within 
10, OOOk0e o of the predicted value. 



8 Appendix The algorithm described in this paper actually measures the probability which is approximately 



2 



proportional to |(3e| . This will give two estimates for e - one a positive range and the other a negative range. Just 
from this measurement, it is not possible to tell the sign of e . However, the sign of e may be easily deduced by 
repeating the algorithm with a slightly perturbed value of the threshold (0. and observing the changes in the estimates 
of e . Denote the subroutine that returns this estimate for e (with the proper sign), given £ Q , & |i by 

e-est e Q , 9) . 

The notation used in this appendix is the same as that used in the paper and discussed in section 2(b), i.e. (J. is the 
threshold, e is the fractional difference in the number of elements with values above and below |i, is the precision 
with which e needs to be estimated. 

(i) Estimating the number of elements above a general threshold The problem discussed previously 

in the main paper assumes that both (0. and £q are given. In case, just the threshold (J. is given but not £q and , the 

algorithm described in this paper could still be used by keeping fixed as a small quantity (e.g. 0.1) and initially 
choosing a high value of Eq (e.g. 0.1) and reducing it by a constant factor in each repetition. The value of the estimate 

that the subroutine e-est (|I, £g, 0) returns will be close to zero until the hypothesized £q becomes of the same 

order as the actual e . The following paragraph gives a more formal description of e-est £ m;M ) • 

Denote this subroutine that returns this estimate for e , given only e & u. by e-est (u, 8 . ) ; this calls the 

° J mm ~ J v ~ min 

previously defined subroutine e-est ((0., £„, 0) a logarithmic number of times. 



{ 




// Note that this subroutine has just two parameters, 
// unlike the one discussed earlier, which has three. 



= 0.1, e Q = 0.1 
{ 



// It is assumed that both e & e 



mm 



are greater than 0.1. 



est = e-est £q, 0) 

if (\est\> 0.2e Q ) break; 

else e n = \est\ X 0.5 



} 




return (est) 



} 



(ii) Estimating the median The median is estimated by using the routine of (i), i.e. e-est E • ) , as a sub- 
routine in a binary search algorithm. It requires a logarithmic number of repetitions of e-est E • ) • 

median (min, max, e . ) 

mm' 

{ 

lower = min; 

upper = max; 

while ( upper - lower > 1 ) 

{ 

upper + lower 

est = e-est (ll,e . ) 

if (est > 0) upper = p, 

else lower = p, 

} 

return (\l) 
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